%load_ext autoreload
%autoreload 2
import math
import additive2.utility as u
import additive2.clusters as cl
import matplotlib.pyplot as plt
import joblib # type: ignore
import pandas as pd
import numpy as np
import seaborn as sns
u.custom_matplotlib_style()
from IPython.display import display # type: ignore
import additive2.notebooks.pore_preformance as nb
The autoreload extension is already loaded. To reload it, use: %reload_ext autoreload
Now let's see how the detected pores differ in terms of pixels. That is if a pixel is detected as pore pixel in the high-res, is it also detected in low res?
# break point
pores_mid_res = joblib.load("../data/Additive/pores-mid-res-v01.np")
pores_high_res = joblib.load("../data/Additive/pores-high-res-v01.np")
pores_pred = joblib.load("../data/Additive/pores-pred-v01.np")
mask = joblib.load("../data/Additive/mask-v01.np")
# let's try to put the three in a single plot
fig = plt.figure(figsize=plt.figaspect(0.3))
axes = []
for i in range(1, 4):
axes.append(fig.add_subplot(1, 3, i, projection='3d'))
pore_imgs = [pores_mid_res, pores_pred, pores_high_res]
titles = ['Mid-resolution', 'Improved resolution', 'High-resolution']
for ax, pore_img, title in zip(axes, pore_imgs, titles):
u.draw_pores(pore_img, mask, ax=ax)
ax.set_title(title)
ax.axis('off')
let's find a way to compute pairwise pixel distance
clusters_high_res, indices_high_res = cl.get_clusters(pores_high_res)
clusters_mid_res, indices_mid_res = cl.get_clusters(pores_mid_res)
clusters_pred, indices_predicted = cl.get_clusters(pores_pred)
indices_high_res.shape
(1711, 3)
clusters_mid, pts_mid, clusters_mid_ = cl.get_pores(pores_mid_res)
clusters_mid_res.max(), clusters_pred.max(), clusters_high_res.max()
(1390, 1865, 3013)
clusters_pred, pts_pred, clusters_pred_ = cl.get_pores(pores_pred)
clusters_high, pts_high, clusters_high_ = cl.get_pores(pores_high_res)
assert len(clusters_mid) == len(pts_mid)
assert len(clusters_high) == len(pts_high)
assert len(clusters_pred) == len(pts_pred)
pred_high_map, high_pred_map = cl.find_matching_pores(pts_pred, pts_high)
mid_high_map, high_mid_map = cl.find_matching_pores(pts_mid, pts_high)
# the following should match pixel-wise comparisons above
(pred_high_map>-1).sum() / len(pred_high_map), (high_pred_map>-1).sum() / len(high_pred_map), (mid_high_map>-1).sum() / len(mid_high_map), (high_mid_map>-1).sum() / len(high_mid_map)
(0.7715996578272027, 0.527177089421391, 0.5581151832460733, 0.3115137346580947)
# the following should match pixel-wise comparisons above
def get_precision_recall(low, high):
tp = np.sum(low>-1)
fp = np.sum(low==-1)
fn = np.sum(high==-1)
return tp/(tp+fp), tp/(tp+fn)
precision_recall_df = pd.DataFrame(
[get_precision_recall(pred_high_map, high_pred_map), get_precision_recall(mid_high_map, high_mid_map)],
columns = ['Precision', 'Recall']
)
precision_recall_df.index = ['Improved Resolution', 'Medium Resolution']
precision_recall_df['Accuracy'] = 1
precision_recall_df.round(2)
| Precision | Recall | Accuracy | |
|---|---|---|---|
| Improved Resolution | 0.77 | 0.53 | 1 |
| Medium Resolution | 0.56 | 0.31 | 1 |
Which it does.
Now let's compare the three images based on the numbers of matching pores
cl_pred_high_map = cl.find_matching_clusters(pred_high_map, clusters_pred, clusters_high)
cl_mid_high_map = cl.find_matching_clusters(mid_high_map, clusters_mid, clusters_high)
# total number of matched pores between mid-res and high res
len(cl_mid_high_map)
14
# total number of matched pores between predicted and high res
len(cl_pred_high_map)
# let's say it is 25
18
# total number of pores detected by high-res is:
np.max(clusters_high)
38
let's check how individual pore sizes compare to each other
# comparing mid vs high pores sizes
res = u.get_pore_size_diff(cl_mid_high_map, clusters_mid_, clusters_high_,
name1='mid', name2='high')
display(res)
res['rel. error'].describe()
| mid | high | rel. error | |
|---|---|---|---|
| 0 | 47 | 54 | 0.15 |
| 1 | 36 | 72 | 1.00 |
| 2 | 39 | 51 | 0.31 |
| 3 | 47 | 61 | 0.30 |
| 4 | 56 | 57 | 0.02 |
| 5 | 37 | 37 | 0.00 |
| 6 | 35 | 61 | 0.74 |
| 7 | 45 | 78 | 0.73 |
| 8 | 44 | 68 | 0.55 |
| 9 | 63 | 79 | 0.25 |
| 10 | 25 | 45 | 0.80 |
| 11 | 33 | 76 | 1.30 |
| 12 | 23 | 46 | 1.00 |
| 13 | 57 | 60 | 0.05 |
count 14.000000 mean 0.514286 std 0.419170 min 0.000000 25% 0.175000 50% 0.430000 75% 0.785000 max 1.300000 Name: rel. error, dtype: float64
res = u.get_pore_size_diff(cl_pred_high_map, clusters_pred_,
clusters_high_, name1='pred', name2='high')
display(res)
res['rel. error'].describe()
| pred | high | rel. error | |
|---|---|---|---|
| 0 | 52 | 55 | 0.06 |
| 1 | 56 | 54 | 0.04 |
| 2 | 66 | 72 | 0.09 |
| 3 | 53 | 51 | 0.04 |
| 4 | 63 | 61 | 0.03 |
| 5 | 37 | 43 | 0.16 |
| 6 | 56 | 57 | 0.02 |
| 7 | 39 | 37 | 0.05 |
| 8 | 50 | 61 | 0.22 |
| 9 | 73 | 78 | 0.07 |
| 10 | 67 | 68 | 0.01 |
| 11 | 42 | 61 | 0.45 |
| 12 | 28 | 81 | 1.89 |
| 13 | 75 | 79 | 0.05 |
| 14 | 36 | 45 | 0.25 |
| 15 | 63 | 76 | 0.21 |
| 16 | 38 | 46 | 0.21 |
| 17 | 55 | 60 | 0.09 |
count 18.000000 mean 0.218889 std 0.431903 min 0.010000 25% 0.042500 50% 0.080000 75% 0.210000 max 1.890000 Name: rel. error, dtype: float64
Getting pore features
# pts_mid
pts_ = pts_mid[np.array(clusters_mid) == 0]
pf = u.PoreFeatures(pts_)
mid_features, high_features_mid = u.get_feature_pairs(cl_mid_high_map, pts_mid,
clusters_mid,
pts_high, clusters_high)
pred_features, high_features_pred = u.get_feature_pairs(cl_pred_high_map, pts_pred, clusters_pred,
pts_high, clusters_high)
u.ttest_1samp([1, 2, 3, 4, 5, 6], 0)
Ttest_1sampResult(statistic=4.58257569495584, pvalue=0.00593354451759226)
stats_ = pd.DataFrame(
(high_features_pred.values - pred_features.values),
columns=pred_features.columns
).assign(tmp=1).groupby('tmp').agg(['mean', 'count', 'std',
u.ttest_pvalue]).T.reset_index()
stats = pd.pivot_table(stats_, index='level_0', columns='level_1',
values=1).drop('cluster')
ci95_hi = []
ci95_lo = []
for i in stats.index:
m, c, s, *_ = stats[['mean', 'count', 'std']].loc[i]
ci95_hi.append(m + 1.96*s/math.sqrt(c))
ci95_lo.append(m - 1.96*s/math.sqrt(c))
stats['ci95_lo'] = ci95_lo
stats['ci95_hi'] = ci95_hi
display(stats.round(2).sort_index())
| level_1 | count | mean | std | ttest_pvalue | ci95_lo | ci95_hi |
|---|---|---|---|---|---|---|
| level_0 | ||||||
| Convex Hull Area | 18.0 | 8.58 | 18.07 | 0.06 | 0.23 | 16.92 |
| Convex Hull Volume | 18.0 | 6.38 | 17.03 | 0.13 | -1.49 | 14.25 |
| Ellipsoid Area | 18.0 | 8.47 | 25.25 | 0.17 | -3.20 | 20.14 |
| Ellipsoid Elongation | 18.0 | -0.05 | 0.84 | 0.82 | -0.44 | 0.34 |
| Ellipsoid Volume | 18.0 | 10.46 | 35.26 | 0.23 | -5.83 | 26.75 |
| Pore Elongation | 18.0 | 0.03 | 0.94 | 0.89 | -0.40 | 0.47 |
| Pore Sphericity | 18.0 | 0.00 | 0.04 | 0.59 | -0.01 | 0.02 |
| Pore Volume | 18.0 | 12.72 | 23.33 | 0.03 | 1.94 | 23.50 |
stats_ = pd.DataFrame(
(high_features_mid.values - mid_features.values),
columns=mid_features.columns
).assign(tmp=1).groupby('tmp').agg(['mean', 'count', 'std',
u.ttest_pvalue]).T.reset_index()
stats = pd.pivot_table(stats_, index='level_0', columns='level_1', values=1).drop('cluster')
ci95_hi = []
ci95_lo = []
for i in stats.index:
m, c, s, *_ = stats[['mean', 'count', 'std']].loc[i]
ci95_hi.append(m + 1.96*s/math.sqrt(c))
ci95_lo.append(m - 1.96*s/math.sqrt(c))
stats['ci95_lo'] = ci95_lo
stats['ci95_hi'] = ci95_hi
display(stats.round(2).sort_index())
| level_1 | count | mean | std | ttest_pvalue | ci95_lo | ci95_hi |
|---|---|---|---|---|---|---|
| level_0 | ||||||
| Convex Hull Area | 14.0 | 32.91 | 39.55 | 0.01 | 12.20 | 53.63 |
| Convex Hull Volume | 14.0 | 14.86 | 17.65 | 0.01 | 5.61 | 24.10 |
| Ellipsoid Area | 14.0 | 37.26 | 60.40 | 0.04 | 5.62 | 68.90 |
| Ellipsoid Elongation | 14.0 | 1.02 | 1.67 | 0.04 | 0.15 | 1.90 |
| Ellipsoid Volume | 14.0 | 21.08 | 34.09 | 0.04 | 3.22 | 38.93 |
| Pore Elongation | 14.0 | 1.36 | 1.70 | 0.01 | 0.47 | 2.25 |
| Pore Sphericity | 14.0 | -0.02 | 0.03 | 0.01 | -0.04 | -0.01 |
| Pore Volume | 14.0 | 38.14 | 32.31 | 0.00 | 21.22 | 55.07 |
# todo write about the model
import additive2.plotEllipsoid as et
plt.figure(figsize=(15, 15))
fig = plt.figure(figsize=(24, 8))
axes = []
for i in range(1, 4):
axes.append(fig.add_subplot(1, 3, i, projection='3d'))
ax = axes[0]
pf = u.PoreFeatures(pts_mid[np.array(clusters_mid)==27])
ellipsoid_tool = et.EllipsoidTool()
u.draw_voxels(pf.pts_in, ax=ax, alpha=.1)
ellipsoid_tool.plotEllipsoid(pf.ellipsoid, plotAxes=True, ax=ax)
ax = axes[2]
pf = u.PoreFeatures(pts_high[np.array(clusters_high)==31])
u.draw_voxels(pf.pts_in, ax=ax, alpha=.1)
ellipsoid_tool = et.EllipsoidTool()
ellipsoid_tool.plotEllipsoid(pf.ellipsoid, plotAxes=True, ax=ax)
ax = axes[1]
pf = u.PoreFeatures(pts_pred[np.array(clusters_pred)==28])
ellipsoid_tool = et.EllipsoidTool()
u.draw_voxels(pf.pts_in, ax=ax, alpha=.1)
ellipsoid_tool.plotEllipsoid(pf.ellipsoid, plotAxes=True, ax=ax)
for ax in axes:
ax.set_xlim([392, 406])
ax.set_ylim([418, 432])
ax.set_zlim([222, 233])
titles = ['Mid-resolution', 'Improved resolution', 'High-resolution']
for ax, title in zip(axes, titles):
ax.set_title(title)
<Figure size 1080x1080 with 0 Axes>
predHighPoreMap_ = pd.DataFrame([list(x) for x in cl_pred_high_map], columns=['pred', 'high'])
midHighPoreMap = pd.DataFrame([list(x) for x in cl_mid_high_map], columns=['mid', 'high'])
predHighPoreMap = pd.concat([
predHighPoreMap_,
pd.DataFrame(set(range(len(clusters_high_))) - set(predHighPoreMap_['high']), columns=['high']),
pd.DataFrame(set(range(len(clusters_pred_))) - set(predHighPoreMap_['pred']), columns=['pred']),
])
common_pores = predHighPoreMap.merge(midHighPoreMap, how='left')
common_pores = pd.concat([
common_pores,
pd.DataFrame(set(range(len(clusters_mid_))) - set(common_pores['mid']), columns=['mid']),
])
common_pores
| pred | high | mid | |
|---|---|---|---|
| 0 | 0.0 | 0.0 | NaN |
| 1 | 2.0 | 2.0 | 1.0 |
| 2 | 3.0 | 9.0 | 7.0 |
| 3 | 5.0 | 10.0 | 8.0 |
| 4 | 6.0 | 12.0 | 10.0 |
| ... | ... | ... | ... |
| 16 | NaN | NaN | 26.0 |
| 17 | NaN | NaN | 28.0 |
| 18 | NaN | NaN | 31.0 |
| 19 | NaN | NaN | 32.0 |
| 20 | NaN | NaN | 33.0 |
74 rows × 3 columns
import additive2.plotEllipsoid as et
plt.figure(figsize=(15, 15))
fig = plt.figure(figsize=(24, 8))
axes = []
for i in range(1, 4):
axes.append(fig.add_subplot(1, 3, i, projection='3d'))
pred_cluster_num, high_cluster_num, mid_cluster_num = common_pores.loc[14]
ax = axes[0]
pf = u.PoreFeatures(pts_mid[clusters_mid_[mid_cluster_num]])
ellipsoid_tool = et.EllipsoidTool()
u.draw_voxels(pf.pts_in, ax=ax, alpha=.1)
ellipsoid_tool.plotEllipsoid(pf.ellipsoid, plotAxes=True, ax=ax)
ax = axes[2]
pf = u.PoreFeatures(pts_high[clusters_high_[high_cluster_num]])
u.draw_voxels(pf.pts_in, ax=ax, alpha=.1)
ellipsoid_tool = et.EllipsoidTool()
ellipsoid_tool.plotEllipsoid(pf.ellipsoid, plotAxes=True, ax=ax)
ax = axes[1]
pf = u.PoreFeatures(pts_pred[clusters_pred_[pred_cluster_num]])
ellipsoid_tool = et.EllipsoidTool()
u.draw_voxels(pf.pts_in, ax=ax, alpha=.1)
ellipsoid_tool.plotEllipsoid(pf.ellipsoid, plotAxes=True, ax=ax)
for ax in axes:
ax.set_xlim([392, 406])
ax.set_ylim([418, 432])
ax.set_zlim([222, 233])
titles = ['Mid-resolution', 'Improved resolution', 'High-resolution']
for ax, title in zip(axes, titles):
ax.set_title(title)
<Figure size 1080x1080 with 0 Axes>
import additive2.plotEllipsoid as et
dTmp = common_pores.dropna()
for j in dTmp.index:
fig = plt.figure(figsize=(24, 8))
axes = []
for i in range(1, 4):
axes.append(fig.add_subplot(1, 3, i, projection='3d'))
pred_cluster_num, high_cluster_num, mid_cluster_num = dTmp.loc[j]
ax = axes[0]
pf_mid = u.PoreFeatures(pts_mid[np.array(clusters_mid)==mid_cluster_num])
ellipsoid_tool = et.EllipsoidTool()
u.draw_voxels(pf_mid.pts_in, ax=ax, alpha=.1)
ellipsoid_tool.plotEllipsoid(pf_mid.ellipsoid, plotAxes=True, ax=ax)
ax = axes[2]
pf_high = u.PoreFeatures(pts_high[np.array(clusters_high)==high_cluster_num])
u.draw_voxels(pf_high.pts_in, ax=ax, alpha=.1)
ellipsoid_tool = et.EllipsoidTool()
ellipsoid_tool.plotEllipsoid(pf_high.ellipsoid, plotAxes=True, ax=ax)
ax = axes[1]
pf_pred = u.PoreFeatures(pts_pred[np.array(clusters_pred)==pred_cluster_num])
ellipsoid_tool = et.EllipsoidTool()
u.draw_voxels(pf_pred.pts_in, ax=ax, alpha=.2)
ellipsoid_tool.plotEllipsoid(pf_pred.ellipsoid, plotAxes=True, ax=ax)
pf_pts = np.concatenate([pf_pred.pts, pf_high.pts, pf_mid.pts])
mn = pf_pts.min(axis=0) - 2
mx = pf_pts.max(axis=0) + 2
for ax in axes:
ax.set_xlim([mn[0], mx[0]])
ax.set_ylim([mn[1], mx[1]])
ax.set_zlim([mn[2], mx[2]])
ax.set_xlabel('x')
ax.set_ylabel('y')
titles = ['Mid-resolution', 'Improved resolution', 'High-resolution']
for ax, title in zip(axes, titles):
ax.set_title(title)
plt.show()
fig, axes = plt.subplots(1, 3, figsize=(15, 5))
axes = axes.reshape(-1)
axes[0].hist([len(x) for x in clusters_mid_], bins=5)
axes[0].set_title('Mid-resolution')
axes[1].hist([len(x) for x in clusters_high_], bins=5)
axes[1].set_title('high-resolution')
axes[2].hist([len(x) for x in clusters_pred_], bins=5);
axes[2].set_title('super-resolution')
for ax in axes:
ax.set_xlim(10, 110)
midPoreFeatures = pd.DataFrame([u.PoreFeatures(pts_mid[x]).all_features for x in clusters_mid_])\
.reset_index().rename(columns={'index': 'pore'})
# midPoreFeatures.columns = ['Mid ' + x for x in midPoreFeatures.columns]
highPoreFeatures = pd.DataFrame([u.PoreFeatures(pts_high[x]).all_features for x in clusters_high_])\
.reset_index().rename(columns={'index': 'pore'})
# highPoreFeatures.columns = ['High ' + x for x in highPoreFeatures.columns]
predPoreFeatures = pd.DataFrame([u.PoreFeatures(pts_pred[x]).all_features for x in clusters_pred_])\
.reset_index().rename(columns={'index': 'pore'})
# predPoreFeatures.columns = ['Super ' + x for x in predPoreFeatures.columns]
len(clusters_mid_), len(clusters_high_), len(clusters_pred_)
(35, 39, 32)
commonPoresUnstacked = common_pores.unstack().rename('pore').reset_index()\
.rename(columns={'level_0': 'type', 'level_1': 'main_pore'})
commonPoresUnstacked.head(4)
| type | main_pore | pore | |
|---|---|---|---|
| 0 | pred | 0 | 0.0 |
| 1 | pred | 1 | 2.0 |
| 2 | pred | 2 | 3.0 |
| 3 | pred | 3 | 5.0 |
commonPoresUnstacked.groupby('type').max().sum()
main_pore 156.0 pore 103.0 dtype: float64
allPoreFeatures = pd.concat([
midPoreFeatures.assign(type='mid'),
highPoreFeatures.assign(type='high'),
predPoreFeatures.assign(type='pred')
])
allPoreFeatures.groupby('type')['pore'].max()
type high 38 mid 34 pred 31 Name: pore, dtype: int64
allSharedPoresFeatures = allPoreFeatures.merge(commonPoresUnstacked, on=['type', 'pore'])
allSharedPoresFeatures
| pore | Convex Hull Area | Pore Volume | Convex Hull Volume | Ellipsoid Volume | Ellipsoid Area | Ellipsoid Elongation | Pore Sphericity | Pore Elongation | type | main_pore | |
|---|---|---|---|---|---|---|---|---|---|---|---|
| 0 | 0 | 82.498975 | 62 | 30.5 | 71.452453 | 136.199815 | 4.375653 | 0.572228 | 5.336956 | mid | 0 |
| 1 | 1 | 135.113122 | 128 | 54.0 | 110.460570 | 204.864846 | 6.230679 | 0.511347 | 7.479158 | mid | 1 |
| 2 | 2 | 81.467619 | 70 | 30.0 | 71.158961 | 135.749488 | 4.368205 | 0.573121 | 5.312026 | mid | 1 |
| 3 | 3 | 77.659737 | 66 | 28.0 | 64.064661 | 124.088692 | 4.116301 | 0.574196 | 5.282254 | mid | 2 |
| 4 | 4 | 97.052389 | 74 | 36.5 | 99.735116 | 185.903807 | 5.711054 | 0.548285 | 6.067092 | mid | 3 |
| ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... | ... |
| 101 | 27 | 99.680190 | 82 | 37.0 | 91.263904 | 172.922892 | 5.489174 | 0.538695 | 6.396917 | pred | 52 |
| 102 | 28 | 107.912687 | 104 | 41.5 | 89.201645 | 168.122077 | 5.280513 | 0.537168 | 6.451623 | pred | 14 |
| 103 | 29 | 208.246167 | 176 | 86.0 | 173.643681 | 317.571353 | 9.391903 | 0.452452 | 10.796470 | pred | 15 |
| 104 | 30 | 119.685728 | 110 | 47.0 | 97.843557 | 182.623290 | 5.625377 | 0.526228 | 6.862447 | pred | 16 |
| 105 | 31 | 145.714777 | 144 | 58.5 | 126.855739 | 233.815149 | 7.023366 | 0.500132 | 7.993679 | pred | 17 |
106 rows × 11 columns
plt.subplots(figsize=(16, 8))
sns.barplot(data=allSharedPoresFeatures,
x='main_pore', y='Pore Volume', hue='type', hue_order=['mid', 'pred', 'high'])
<AxesSubplot:xlabel='main_pore', ylabel='Pore Volume'>
plt.figure(figsize=(18, 12))
plt.subplot(212)
colName = 'Pore Volume'
colNameNew = colName + r" ($\mu m^3$)"
dTmp = nb.GetMatchedOrOnMatched(allSharedPoresFeatures, common_pores, 'high', 'pred')
dTmp[colNameNew] = dTmp[colName]
sns.scatterplot(data=dTmp.drop([14]), x='pore', y=colNameNew, hue='Is Detected')
plt.title("(b) Pores detected in improved-resolution data")
plt.subplot(211)
dTmp = nb.GetMatchedOrOnMatched(allSharedPoresFeatures, common_pores, 'high', 'mid')
dTmp[colNameNew] = dTmp[colName]
sns.scatterplot(data=dTmp, x='pore', y=colNameNew, hue='Is Detected')
plt.title("(a) Pores detected in mid-resolution data")
/home/bzr0014/git/additive2/additive2/notebooks/pore_preformance.py:13: SettingWithCopyWarning: A value is trying to be set on a copy of a slice from a DataFrame. Try using .loc[row_indexer,col_indexer] = value instead See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy tmp['has_match'] = ~tmp[predCol].isna() /home/bzr0014/git/additive2/additive2/notebooks/pore_preformance.py:13: SettingWithCopyWarning: A value is trying to be set on a copy of a slice from a DataFrame. Try using .loc[row_indexer,col_indexer] = value instead See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy tmp['has_match'] = ~tmp[predCol].isna()
Text(0.5, 1.0, '(a) Pores detected in mid-resolution data')
filt = dTmp['Is Detected']
from scipy.stats import ttest_ind
ttest_ind(dTmp.loc[filt, 'Pore Volume'], dTmp.loc[~filt, 'Pore Volume'])
Ttest_indResult(statistic=5.7529138721912005, pvalue=1.3569643216397724e-06)
from scipy.spatial.distance import pdist, cdist
a = []
b = []
c = []
for _, (i, j) in common_pores[['high', 'pred']].dropna().iterrows():
tmp1 = pts_high[clusters_high_[i]]
tmp2 = pts_pred[clusters_pred_[j]]
tmp = cdist(tmp1, tmp2) == 0
a.append(tmp.sum())
b.append(tmp1.shape[0])
c.append(tmp2.shape[0])
(
np.array(a).sum() / np.array(b).sum(),
np.array(a).sum() / np.array(c).sum(),
)
(0.8313364055299539, 0.9504741833508957)
a = []
b = []
c = []
for _, (i, j) in common_pores[['high', 'mid']].dropna().iterrows():
tmp1 = pts_high[clusters_high_[i]]
tmp2 = pts_mid[clusters_mid_[j]]
tmp = cdist(tmp1, tmp2) == 0
a.append(tmp.sum())
b.append(tmp1.shape[0])
c.append(tmp2.shape[0])
(
np.array(a).sum() / np.array(b).sum(),
np.array(a).sum() / np.array(c).sum(),
)
(0.6307692307692307, 0.9080068143100511)
def elongation(r, l):
r = np.reshape(r, (-1, 1))
l = np.reshape(l, (1, -1))
a = 2 * np.pi * r * l
v = np.pi * r * r * l
return a**3/v**2
x = np.arange(1, 10)#.round(2)
res = elongation(x, x)
fig, ax = plt.subplots(figsize=(10, 10))
x_, y_ = np.meshgrid(x, x)
ax.imshow(np.log(res+1), origin='lower')
ax.set_xticks(x);
ax.set_xticklabels(x);
ax.set_yticklabels(x);
ax.set_xlabel('L')
ax.set_ylabel('R')
/home/bzr0014/miniconda3/envs/additive2/lib/python3.6/site-packages/ipykernel_launcher.py:6: UserWarning: FixedFormatter should only be used together with FixedLocator
Text(0, 0.5, 'R')
elongation(10, 1), elongation(1, 10), elongation(1, 1), elongation(10, 10)
(array([[2.51327412]]), array([[251.32741229]]), array([[25.13274123]]), array([[25.13274123]]))
plt.plot(res[:, 5])
[<matplotlib.lines.Line2D at 0x7fa3d1c2d7f0>]